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The  Automatic-Flow  (AFLOW)  standard  for  the  high-throughput  construction  of  materials  science  elec¬ 
tronic  structure  databases  is  described.  Electronic  structure  calculations  of  solid  state  materials  depend 
on  a  large  number  of  parameters  which  must  be  understood  by  researchers,  and  must  be  reported  by 
originators  to  ensure  reproducibility  and  enable  collaborative  database  expansion.  We  therefore  describe 
standard  parameter  values  for  k-point  grid  density,  basis  set  plane  wave  kinetic  energy  cut-off, 
exchange-correlation  functionals,  pseudopotentials,  DFT+U  parameters,  and  convergence  criteria  used 
in  AFLOW  calculations. 
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1.  Introduction 

The  emergence  of  computational  materials  science  over  the  last 
two  decades  has  been  inextricably  linked  to  the  development  of 
complex  quantum-mechanical  codes  that  enable  accurate  evalua¬ 
tion  of  the  electronic  and  thermodynamic  properties  of  a  wide 
range  of  materials.  The  continued  advancement  of  this  field  entails 
the  construction  of  large  open  databases  of  materials  properties 
that  can  be  easily  reproduced  and  extended.  One  obstacle  to  the 
reproducibility  of  the  data  is  the  unavoidable  complexity  of  the 
codes  used  to  obtain  it.  Published  data  usually  includes  basic 
information  about  the  underlying  calculations  that  allows  rough 
reproduction.  However,  exact  duplication  depends  on  many  details, 
that  are  seldom  reported,  and  is  therefore  difficult  to  achieve. 

These  difficulties  might  limit  the  utility  of  the  databases  cur¬ 
rently  being  created  by  high-throughput  frameworks,  such  as 
AFLOW  [1-3]  and  the  Materials  Project  [4,5].  For  maximal  impact, 
the  data  stored  in  these  repositories  must  be  generated  and  repre¬ 
sented  in  a  consistent  and  robust  manner,  and  shared  through 
standardized  calculation  and  communication  protocols.  Following 
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these  guidelines  would  promote  optimal  use  of  the  results  gener¬ 
ated  by  the  entire  community. 

The  AFLOW  (Automatic  FLOW)  code  is  a  framework  for 
high-throughput  computational  materials  discovery  [1-3,6],  using 
separate  DFT  packages  to  calculate  electronic  structure  and  opti¬ 
mize  the  atomic  geometry.  The  AFLOW  framework  works  with 
the  VASP  [7-10]  DFT  package,  and  integration  with  the  Quantum 
ESPRESSO  software  [11]  is  currently  in  progress.  The  AFLOW 
framework  includes  preprocessing  functions  for  generating  input 
files  for  the  DFT  package;  obtaining  the  initial  geometric  structures 
by  extracting  the  relevant  data  from  crystallographic  information 
files  or  by  generating  them  using  inbuilt  prototype  databases, 
and  then  transforming  them  into  standard  forms  which  are  easiest 
to  calculate.  It  then  runs  and  monitors  the  DFT  calculations  auto¬ 
matically,  detecting  and  responding  to  calculation  failures, 
whether  they  are  due  to  insufficient  hardware  resources  or  to  run¬ 
time  errors  of  the  DFT  calculation  itself.  Finally,  AFLOW  contains 
postprocessing  routines  to  extract  specific  properties  from  the 
results  of  one  or  more  of  the  DFT  calculations,  such  as  the  band 
structure  or  thermal  properties  [12], 

The  AFLOWLIB  repository  [2,3,6]  was  built  according  to  these 
principles  of  consistency  and  reproducibility,  and  the  data  it  con¬ 
tains  can  be  easily  accessed  through  a  REpresentational  State 
Transfer-Application  Programming  Interface  (REST-AP1)  [3[.  In  this 
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paper  we  present  a  detailed  description  of  the  AFLOW  standard  for 
high-throughput  (HT)  materials  science  calculations  by  which  the 
data  in  this  repository  was  created. 


2.  AFLOW  calculation  types 


the  Automatic  Gibbs  Library  (AGL)  [12],  and  the  Automatic 
Phonon  Library  (APL)  [1  j,  which  are  methods  that  have  been  imple¬ 
mented  within  the  AFLOW  framework.  In  the  following,  we 
describe  the  parameter  sets  used  to  address  the  particular  chal¬ 
lenges  of  the  calculations  included  in  each  AFLOW  repository. 


The  AFLOWL1B  consortium  [2]  repository  is  divided  into  data¬ 
bases  containing  calculated  properties  of  over  625,000  materials: 
the  Binary  Alloy  Project,  the  Electronic  Structure  database,  the 
Heusler  database,  and  the  Elements  database.  These  are  freely 
accessible  online  via  the  AFLOWLIB  website  [6  ,  as  well  as  through 
the  API  [3],  The  Electronic  Structure  database  consists  of  entries 
found  in  the  Inorganic  Crystal  Structures  Database,  ICSD  [13,14], 
and  will  thus  be  referred  to  as  “ICSD”  throughout  this  publication. 
The  Heusler  database  consists  of  ternary  compounds,  primarily 
based  on  the  Heusler  structure  but  with  other  structure  types 
now  being  added. 

The  high-throughput  construction  of  these  materials  databases 
relies  on  a  pre-defined  set  of  standard  calculation  types.  These  are 
designed  to  accommodate  the  interest  in  various  properties 
of  a  given  material  (e.g.  the  ground  state  ionic  configuration, 
thermodynamic  quantities,  electronic  and  magnetic  properties), 
the  program  flow  of  the  HT  framework  that  envelopes  the  DFT 
portions  of  the  calculations,  as  well  as  the  practical  need  for 
computational  robustness.  The  AFLOW  standard  thus  deals  with 
the  parameters  involved  in  the  following  calculation  types: 

i.  RELAX.  Geometry  optimizations  using  algorithms  imple¬ 
mented  within  the  DFT  package.  This  calculation  type  is  con¬ 
cerned  with  obtaining  the  ionic  configuration  and  cell  shape 
and  volume  that  correspond  to  a  minimum  in  the  total 
energy.  It  consists  of  two  sequential  relaxation  steps.  The 
starting  point  for  the  first  step,  RELAXl,  can  be  an  entry 
taken  from  an  external  source,  such  as  a  library  of  alloy  pro¬ 
totypes  [15,16],  the  ICSD  database,  or  the  Pauling  File  [17], 
These  initial  entries  are  preprocessed  by  AFLOW,  and  cast 
into  a  unit  cell  that  is  most  convenient  for  calculation,  usu¬ 
ally  the  standard  primitive  cell,  in  the  format  appropriate 
for  the  DFT  package  in  use.  The  second  step,  RELAX2,  uses 
the  final  ionic  positions  from  the  first  step  as  its  starting 
point,  and  serves  as  a  type  of  annealing  step.  This  is  used 
for  jumping  out  of  possible  local  minima  resulting  from 
wavefunction  artifacts. 

ii.  STATIC.  A  single-point  energy  calculation.  The  starting 
point  is  the  set  of  final  ionic  positions,  as  produced  by  the 
RELAX 2  step.  The  outcome  of  this  calculation  is  used  in  the 
determination  of  most  of  the  thermodynamic  and  electronic 
properties  included  in  the  various  AFLOW  databases.  It 
therefore  applies  a  more  demanding  set  of  parameters  than 
those  used  on  the  RELAX  set  of  runs. 

iii.  BANDS.  Electronic  band  structure  generation.  The  converged 
STATIC  charge  density  and  ionic  positions  are  used  as  the 
starting  points,  and  the  wavefunctions  are  reoptimized 
along  standardized  high  symmetry  lines  connecting  special 
fc-points  in  the  irreducible  Brillouin  zone  (IBZ)  [18], 

These  calculation  types  are  performed  in  the  order  shown  above 
(i.e.  RELAXl  RELAX2  ->  STATIC  BANDS)  on  all  materials 
found  in  the  Elements,  ICSD,  and  Heusler  databases.  Those  found 
in  the  Binary  Alloy  database  contain  data  produced  only  by  the 
two  RELAX  calculations.  Sets  of  these  calculation  types  can  be  com¬ 
bined  to  describe  more  complex  phenomena  than  can  be  obtained 
from  a  single  calculation.  For  example,  sets  of  RELAX  and  STATIC 
calculations  for  different  cell  volumes  and/or  atomic  configura¬ 
tions  are  used  to  calculate  thermal  and  mechanical  properties  by 


3.  The  AFLOW  standard  parameter  set 

The  standard  parameters  described  in  this  work  are  classified 
according  to  the  wide  variety  of  tasks  that  a  typical  solid  state 
DFT  calculation  involves:  Brillouin  zone  sampling,  Fourier  transform 
meshes,  basis  sets,  potentials,  self-interaction  error  (SIE) 
corrections,  electron  spin,  algorithms  guiding  SCF  convergence 
and  ionic  relaxation,  and  output  options. 

Due  to  the  intrinsic  complexity  of  the  DFT  codes  it  is  impractical 
to  specify  the  full  set  of  DFT  calculation  parameters  within  an  HT 
framework.  Therefore,  the  AFLOW  standard  adopts  many,  but  not 
all,  of  the  internal  defaults  set  by  the  DFT  software  package.  This 
is  most  notable  in  the  description  of  the  Fourier  transform  meshes, 
which  rely  on  a  discretization  scheme  that  depends  on  the  applied 
basis  and  crystal  geometry  for  its  specification.  Those  internal 
default  settings  are  cast  aside  when  error  corrections  of  failed 
DFT  runs,  an  integral  part  of  AFLOW’s  functionality,  take  place. 
The  settings  described  in  this  work  are  nevertheless  prescribed 
as  fully  as  is  practicable,  in  the  interest  of  providing  as  much  infor¬ 
mation  as  possible  to  anyone  interested  in  reproducing  or  building 
on  our  results. 


3.1.  k-point  sampling 

Two  approaches  are  used  when  sampling  the  IBZ:  the  first  con¬ 
sists  of  uniformly  distributing  a  large  number  of  fc-points  in  the  IBZ, 
while  the  second  relies  on  the  construction  of  paths  connecting 
high  symmetry  (special)  fc-points  in  the  IBZ.  Within  AFLOW,  the 
second  sampling  method  corresponds  to  the  BANDS  calculation 
type,  whereas  the  other  calculation  types  (non-BANDS)  are  per¬ 
formed  using  the  first  sampling  method. 

Sampling  in  non-BANDS  calculations  is  obtained  by  defining  and 
setting  NKPPRA,  the  number  of  fc-points  per  reciprocal  atom.  This 
quantity  determines  the  total  number  of  fc-points  in  the  IBZ,  taking 
into  account  the  fc-points  density  along  each  reciprocal  lattice  vec¬ 
tor  as  well  as  the  number  of  atoms  in  the  simulation  cell,  via  the 
relation: 


^kppra  ^  tnin 


n* 


xNa 


(1) 


Na  is  the  number  of  atoms  in  the  cell,  and  the  IV,  factors  correspond 
to  the  number  of  sampling  points  along  each  reciprocal  lattice  vec¬ 
tor,  bj,  respectively.  These  factors  define  the  grid  resolution, 
<5k,j|f>j||/JVj,  which  is  made  as  uniform  as  possible  under  the  con¬ 
straint  of  Eq.  (1).  The  fc-point  meshes  are  then  generated  within 
the  Monkhorst-Pack  scheme  [19],  unless  the  material  belongs  to 
the  hP,  or  hR  Bravais  lattices,  in  which  case  the  hexagonal  symmetry 
is  preserved  by  centering  the  mesh  at  the  r -point. 

Default  Nkppra  values  depend  on  the  calculation  type  and  the 
database.  The  NKPPRA  values  used  for  the  entries  in  the  Elements 
database  are  material  specific  and  set  manually  due  to  convergence 
of  the  total  energy  calculation.  The  defaults  applied  to  the  RELAX 
and  STATIC  calculations  are  summarized  in  Table  1.  These  defaults 
ensure  proper  convergence  of  the  calculations.  They  may  be  too 
stringent  for  some  cases  but  enable  reliable  application  within 
the  HT  framework,  thus  presenting  a  practicable  balance  between 
accuracy  and  calculation  cost. 
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Table  1 

Default  Nkppm  values  used  in  non-BAUDS  calculations. 


Database 

STATIC 

RELAX 

Binary  Alloy 

N.A. 

6000 

Heusler 

10,000 

6000 

ICSD 

10,000 

8000 

For  BANDS  calculations  AFLOW  generates  Brillouin  zone  inte¬ 
gration  paths  in  the  manner  described  in  a  previous  publication 
[18].  The  k-point  sampling  density  is  the  line  density  of  k-points 
along  each  of  the  straight-line  segments  of  the  path  in  the  IBZ. 
The  default  setting  of  AFLOW  is  128  k-points  along  each  segment 
connecting  high-symmetry  k-points  in  the  IBZ  for  single  element 
structures,  and  20  k-points  for  compounds. 

The  occupancies  at  the  Fermi  edge  in  all  non-RELAX  type  runs 
are  handled  via  the  tetrahedron  method  with  Blochl  corrections 
[20].  This  involves  the  NKPPRA  parameter,  as  described  above,  in 
RELAX  type  calculations,  where  the  determination  of  accurate 
forces  is  important,  some  type  of  smearing  must  be  performed. 
In  cases  where  the  material  is  assumed  to  be  a  metal,  the 
Methfessel-Paxton  approach  [21]  is  adopted,  with  a  smearing 
width  of  0.10  eV.  Gaussian  smearing  is  used  in  all  other  types  of 
materials,  with  a  smearing  width  of  0.05  eV. 


3.2.  Potentials  and  basis  set 

The  interactions  involving  the  valence  electron  shells  are  han¬ 
dled  with  the  potentials  provided  with  the  DFT  software  package. 
In  VASP,  these  include  Ultra-Soft  Pseudopotentials  (USPP)  [22,23] 
and  Projector-Augmented  Wavefunction  (PAW)  potentials 
[24,25],  which  are  constructed  according  to  the  Local  Density 
Approximation  (LDA)  [26,27],  and  the  Generalized  Gradient 
Approximation  (GGA)  PW91  [28,29]  and  PBE  [30,31]  exchange- 
correlation  (XC)  functionals.  The  ICSD,  Binary  Alloy  and  Heusler 
databases  built  according  to  the  AFLOW  standard  use  the  PBE  func¬ 
tional  combined  with  the  PAW  potential  as  the  default.  The  PBE 
functional  is  among  the  best  studied  GGA  functionals  used  in  crys¬ 
talline  systems,  while  the  PAW  potentials  are  preferred  due  to  their 
advantages  over  the  USPP  methodology.  Nevertheless,  defaults 
have  been  defined  for  a  number  of  potential/XC  functional  combi¬ 
nations,  and  in  the  case  of  the  Elements  database,  results  are  avail¬ 
able  for  LDA,  GGA-PW91  and  GGA-PBE  functionals  with  both  USPP 
and  PAW  potentials.  Additionally,  there  are  a  small  number  of 
entries  in  the  ICSD  and  Binary  Alloy  databases  (less  than  1%  of 
the  total)  which  have  been  calculated  with  the  GGA-PW91  func¬ 
tional  using  either  the  USPP  or  PAW  potential.  The  exact  combina¬ 
tion  of  exchange-correlation  functional  and  potential  used  for  a 
specific  entry  in  the  AFLOWLIB  database  can  always  be  determined 
by  querying  the  keyword  df  t_type  using  the  AFLOWLIB  REST-API 
[3]. 

DFT  packages  often  provide  more  than  one  potential  of  each 
type  per  element.  The  AFLOW  standardized  lists  of  PAW  and 
USPP  potentials  are  presented  in  Tables  2  and  3,  respectively.  The 
“Label”  column  in  these  tables  corresponds  to  the  naming  conven¬ 
tion  adopted  by  VASP.  The  checksum  of  each  file  listed  in  the  tables 
is  included  in  the  accompanying  supplement  for  verification 
purposes. 

Each  potential  provided  with  the  VASP  package  has  two  recom¬ 
mended  plane-wave  kinetic  energy  cut-off  (£„,)  values,  the  smaller 
of  which  ensures  the  reliability  of  a  calculation  to  within  a 
well-defined  error.  Additionally,  materials  with  more  than  one  ele¬ 
ment  type  will  have  two  or  more  sets  of  recommended  Ecut  values. 
In  the  AFLOW  standard,  the  applied  £n,t  value  is  the  largest  found 
among  the  recommendations  for  all  species  involved  in  the  calcu¬ 
lation,  increased  by  a  factor  of  1 .4. 


Table  2 

Projector-Augmented  Wavefunction  (PAW)  potentials,  parameterized  for  the  LDA, 
PW91,  and  PBE  functionals,  included  in  the  AFLOW  standard.  The  PAW-PBE 
combination  is  used  as  the  default  for  ICSD,  Binary  Alloy  and  Heusler  databases. 


Element 

Label 

Element 

Label 

Element 

Label 

H 

H 

Se 

Se 

Gdb 

Gd_3 

He 

He 

Br 

Br 

Tb 

Tb_3 

Li 

LLsv 

Kr 

Kr 

Dy 

Dy_3 

Be 

Be_sv 

Rb 

Rb_sv 

Ho 

Ho_3 

B 

B_h 

Sr 

Sr_sv 

Er 

Er_3 

C 

c 

Y 

Y_sv 

Tm 

Tm 

N 

N 

Zr 

Zr_sv 

Yb 

Yb 

O 

O 

Nb 

Nb_sv 

Lu 

Lu 

F 

F 

Mo 

Mo_pv 

Hf 

Hf 

Ne 

Ne 

Tc 

Tc_pv 

Ta 

Ta_pv 

Na 

Na_pv 

Ru 

Ru_pv 

W 

W_pv 

Mg 

Mg_pv 

Rh 

Rh_pv 

Re 

Re_pv 

Al 

Al 

Pd 

Pd_pv 

Os 

Os_pv 

Si 

Si 

Ag 

Ag 

Ir 

Ir 

P 

P 

Cd 

Cd 

Pt 

Pt 

S 

S 

In 

In_d 

Au 

Au 

Cl 

Cl 

Sn 

Sn 

Hg 

Hg 

Ar 

Ar 

Sb 

Sb 

Tl 

TLd 

K 

K_sv 

Te 

Te 

Pb 

Pb_d 

Ca 

Ca_sv 

I 

I 

Bi 

BLd 

Sc 

Sc_sv 

Xe 

Xe 

Po 

Po 

Ti 

TLsv 

Cs 

Cs_sv 

At 

At 

V 

V_sv 

Ba 

Ba_sv 

Rn 

Rn 

Cr 

Cr_pv 

La 

La 

Fr 

Fr 

Mn 

Mn_pv 

Ce 

Ce 

Ra 

Ra 

Fe 

Fe_pv 

Pr 

Pr 

Ac 

Ac 

Co 

Co 

Nd 

Nd 

Th 

Th_s 

Ni 

NLpv 

Pm 

Pm 

Pa 

Pa 

Cu 

Cu_pv 

Smc1 

Sm 

U 

U 

Zn 

Zn 

Smb 

Sm_3 

Np 

Np_s 

Ga 

Ga_h 

Eu 

Eu 

Pu 

Pu_s 

As 

As 

Gd" 

Gd 

a  PBE  potentials  only. 
b  LDA  and  PW91  potentials  only. 


It  is  possible  to  evaluate  the  non-local  parts  of  the  potentials  in 
real  space,  rather  than  in  the  more  computationally  intensive 
reciprocal  space.  This  approach  is  prone  to  aliasing  errors,  and 
requires  the  optimization  of  real-space  projectors  if  these  are  to 
be  avoided.  The  real-space  projection  scheme  is  most  appropriate 
for  larger  systems,  e.g.  surfaces,  and  is  therefore  not  used  in  the 
construction  of  the  databases  found  in  the  AFLOWLIB  repository. 

3.3.  Fourier  transform  meshes 

As  mentioned  previously,  it  is  not  practical  to  describe  the  pre¬ 
cise  default  settings  that  are  applied  by  the  AFLOW  standard  in  the 
specification  of  the  Fourier  transform  meshes.  We  shall  just  note 
that  they  are  defined  in  terms  of  the  grid  spacing  along  each  of 
the  reciprocal  lattice  vectors,  S,.  These  are  obtained  from  the  set 
of  real  space  lattice  vectors,  a,-,  via  [b3b2b3]  =  27t[aia2a3]  .  A  dis¬ 
tance  in  reciprocal  space  is  then  defined  by  d,  =  ||b,]|/nj,  where 
the  set  of  n,  are  the  number  of  grid  points  along  each  reciprocal  lat¬ 
tice  vector,  and  where  the  total  number  of  points  in  the  simulation 
is  Hi  x  n2  x  n3. 

The  VASP  package  relies  primarily  on  the  so-called  dual  grid 
technique,  which  consists  of  two  overlapping  meshes  with  different 
coarseness.  The  least  dense  of  the  two  is  directly  dependent  on  the 
applied  plane-wave  basis,  Ecut,  while  the  second  is  a  finer  mesh 
onto  which  the  charge  density  is  mapped.  The  AFLOW  standard 
relies  on  placing  sufficient  points  in  the  finer  mesh  such  that 
wrap-around  (“aliasing”)  errors  are  avoided.  In  terms  of  the  quan¬ 
tity  d i,  defined  above,  the  finer  grid  is  characterized  by 
d,  R5  0.10  A-1,  while  the  coarse  grid  results  in  d,  ss  0.15  A-1.  These 
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Table  3 

Ultra-Soft  Pseudopotentials  (USPP),  parameterized  for  the  LDA  and  PW91  functionals, 
included  in  the  AFLOW  standard. 


Element 

Label 

Element 

Label 

Element 

Label 

H 

H_soft 

As 

As 

Tb 

Tb_3 

He 

He 

Se 

Se 

Dy 

Dy_3 

Li 

LLpv 

Br 

Br 

Ho 

Ho_3 

Be 

Be 

Kr 

Kr 

Er 

Er_3 

B 

B 

Rb 

Rb_pv 

Tm 

Tm 

C 

C 

Sr 

Sr_pv 

Yb 

Yb 

N 

N 

Y 

Y_pv 

Lu 

Lu 

0 

0 

Zr 

Zr_pv 

Hf 

Hf 

F 

F 

Nb 

Nb_pv 

Ta 

Ta 

Ne 

Ne 

Mo 

Mo_pv 

W 

W 

Na 

Na_pv 

Tc 

Tc 

Re 

Re 

Mg 

Mg_pv 

Ru 

Ru 

Os 

Os 

Al 

Al 

Rh 

Rh 

Ir 

Ir 

Si 

Si 

Pd 

Pd 

Pt 

Pt 

P 

P 

Ag 

Ag 

Au 

Au 

S 

S 

Cd 

Cd 

Hg 

Hg 

Cl 

Cl 

In 

In_d 

Tl 

TLd 

Ar 

Ar 

Sn 

Sn 

Pb 

Pb 

I< 

K_pv 

Sb 

Sb 

Bi 

Bi 

Ca 

Ca_pv 

Te 

Te 

Po 

Po 

Sc 

Sc_pv 

I 

I 

At 

At 

Ti 

Ti_pv 

Xe 

Xe 

Rn 

Rn 

V 

V_pv 

Cs 

Cs_pv 

Fr 

Fr 

Cr 

Cr 

Ba 

Ba_pv 

Ra 

Ra 

Mn 

Mn 

La 

La 

Ac 

Ac 

Fe 

Fe 

Ce 

Ce 

Th 

Th_s 

Co 

Co 

Pr 

Pr 

Pa 

Pa 

Ni 

Ni 

Nd 

Nd 

U 

U 

Cu 

Cu 

Pm 

Pm 

Np 

Np_s 

Zn 

Zn 

Sm 

Sm_3 

Pu 

Pu_s 

Ga 

Ga_d 

Eu 

Eu 

Ge 

Ge 

Gd 

Gd 

two  values  are  approximate,  as  there  is  significant  dispersion  in 
these  quantities  across  the  various  databases. 

3.4.  DFT+U  corrections 

Extended  systems  containing  d  and  /  block  elements  are  often 
poorly  represented  within  DFT  due  to  the  well  known  self  interac¬ 
tion  error  (S1E)  [27].  The  influence  that  the  SIE  has  on  the  energy 
gap  of  insulators  has  long  been  recognized,  and  several  methods 
that  account  for  it  are  available.  These  include  the  GW  approxima¬ 
tion  [32],  the  rotationally  invariant  approach  introduced  by 
Dudarev  [33]  and  Liechtenstein  [34]  (denoted  here  as  DFT+U),  as 
well  as  the  recently  developed  ACBNO  pseudo-hybrid  density  func¬ 
tional  [35], 

The  DFT+U  approach  is  currently  the  best  suited  for 
high-throughput  investigations,  and  is  therefore  included  in  the 
AFLOW  standard  for  the  entire  ICSD  database,  and  is  also  used 
for  certain  entries  in  the  Heusler  database  containing  the  elements 
0,  S,  Se,  and  F.  It  is  not  used  for  the  Binary  Alloy  database.  This 
method  has  a  significant  dependence  on  parameters,  as  each  atom 
is  associated  with  two  numbers,  the  screened  Coulomb  parameter, 
U,  and  the  Stoner  exchange  parameter,  J.  These  are  usually  reported 
as  a  single  factor,  combined  via  Ueff  =  U  -J.  The  set  of  Ue ff  values 
associated  with  the  d  block  elements  [18,36]  are  presented  in 
Table  4,  to  which  the  elements  In  and  Sn  have  been  added. 

A  subset  of  the  /-block  elements  can  be  found  among  the  sys¬ 
tems  included  in  the  AFLOWLIB  consortium  databases.  We  are 
not  aware  of  the  existence  of  a  systematic  search  for  the  best  set 
of  U  and  J  parameters  for  the  elements  Nd,  Sm,  and  Eu,  so  we  have 
relied  on  an  in-house  parameterization  [18]  for  those  entries  in  the 
databases  that  contain  them.  The  values  used  are  reproduced  in 
Table  5.  Note  that  by  construction  the  SIE  correction  must  be 


Table  4 

t/eff  parameters  applied  to  d  orbitals. 


Element 

Uetf 

Refs. 

Element 

f^eff 

Refs. 

Sc 

2.9 

[38] 

W 

2.2 

[39] 

Ti 

4.4 

[40] 

Tc 

2.7 

[39] 

V 

2.7 

[411 

Ru 

3.0 

[39] 

Cr 

3.5 

[42] 

Rh 

3.3 

[39] 

Mn 

4.0 

[42] 

Pd 

3.6 

[391 

Fe 

4.6 

[431 

Ag 

5.8 

[44] 

Co 

5.0 

[411 

Cd 

2.1 

[45] 

Ni 

5.1 

[411 

In 

1.9 

[45] 

Cu 

4.0 

[42] 

Sn 

3.5 

[46] 

Zn 

7.5 

[451 

Ta 

2.0 

[39] 

Ga 

3.9 

[471 

Re 

2.4 

[39] 

Sn 

3.5 

[46] 

Os 

2.6 

[39] 

Nb 

2.1 

[391 

Ir 

2.8 

[39] 

Mo 

2.4 

[391 

Pt 

3.0 

[39] 

Ta 

2.0 

[461 

Au 

4.0 

Table  5 

U  and  J  parameters  applied  to  selected  /-block  elements. 

Element 

u 

j 

Refs. 

Element 

u 

j 

Refs. 

La 

8.1 

0.6 

[48] 

Dy 

5.6 

0.0 

[49] 

Ce 

7.0 

0.7 

[50] 

Tm 

7.0 

1.0 

[51] 

Pr 

6.5 

1.0 

[52] 

Yb 

7.0 

0.67 

[53] 

Nd 

7.2 

1.0 

[18] 

Lu 

4.8 

0.95 

[48] 

Sm 

7.4 

1.0 

[18] 

Th 

5.0 

0.0 

[54] 

Eu 

6.4 

1.0 

[18] 

u 

4.0 

0.0 

[55] 

Gd 

6.7 

0.1 

[56] 

applied  to  a  pre-selected  value  of  the  ^-quantum  number,  and  all 
elements  listed  in  Table  4  correspond  to  i  =  2,  while  those  found 
in  Table  5  correspond  to  l  =  3. 

The  U  and  J  values  listed  have  only  been  applied  to  neutral  sys¬ 
tems,  given  that  there  are  no  entries  in  the  AFLOWLIB  repositories 
that  contain  determinations  of  charged  states.  The  in-house 
parameterization  [18]  applicable  to  Nd,  Sm,  and  Eu  was  performed 
by  fitting  their  4/ levels  to  the  corresponding  experimental  density 
of  states  obtained  from  X-ray  Photoelectron  Spectroscopy-Brems 
strahlung  Isochromat  Spectroscopy  (XPS-BIS)  data  [37],  All  other 
DFT+U  parameters  have  been  taken  from  the  literature,  and  the 
corresponding  citations  are  also  contained  in  Tables  4  and  5. 

3.5.  Spin  polarization 

The  first  of  the  two  RELAX  calculations  is  always  performed  in  a 
collinear  spin-polarized  fashion.  The  initial  magnetic  moments  in 
this  step  are  set  to  the  number  of  atoms  in  the  system,  e.g. 
1.0  pB/atom.  If  the  magnetization  resulting  from  the  RELAX  1  step 
is  found  to  be  below  0.025  pB/ atom,  AFLOW  economizes  computa¬ 
tional  resources  by  turning  spin  polarization  off  in  all  ensuing  cal¬ 
culations.  Spin-orbit  coupling  is  not  used  in  the  current  AFLOW 
standard,  since  it  is  still  too  expensive  to  include  in  a  HT 
framework. 

3.6.  Calculation  methods  and  convergence  criteria 

Two  nested  loops  are  involved  in  the  DFT  calculations  used  by 
AFLOW  in  the  construction  of  the  databases.  The  inner  loop  con¬ 
tains  routines  that  iteratively  optimize  the  electronic  degrees  of 
freedom  (EDOF),  and  features  a  number  of  algorithms  that  are  con¬ 
cerned  with  diagonalizing  the  Kohn-Sham  (KS)  Hamiltonian  at 
each  iteration.  The  outer  loop  performs  adjustments  to  the  system 
geometry  (ionic  degrees  of  freedom,  IDOF)  until  the  forces  acting 
on  the  system  are  minimized. 
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The  convergence  condition  for  each  loop  has  been  defined  in 
terms  of  an  energy  difference,  SE.  If  successive  energies  resulting 
from  the  completion  of  a  loop  are  denoted  as  £,_i  and  £,-,  then  con¬ 
vergence  is  met  when  the  condition  SE  >  £,  -  £,-  i  is  fulfilled.  Note 
that  £,  can  either  be  the  electronic  energy  resulting  from  the  inner 
loop,  or  the  configurational  energy  resulting  from  the  outer  loop. 
The  electronic  convergence  criteria  will  be  denoted  as  SEeieC,  and 
the  ionic  criteria  as  SEion.  The  AFLOW  standard  relies  on 
dEeiec  =  1 CT5  eV  and  <5£fon  =  1CT4  eV  for  entries  in  the  Elements 
database.  All  other  databases  include  calculations  performed  with 
dEeiec  =  1 0~3  eV  and  <5£ion  =  1  CT2  eV. 

Optimizations  of  the  EDOF  depend  on  sets  of  parameters  that 
fall  under  three  general  themes:  initial  guesses,  diagonalization 
methods,  and  charge  mixing.  The  outer  loop  (optimizations  of 
the  1DOF)  is  concerned  with  the  lattice  vectors  and  the  ionic  posi¬ 
tions,  and  is  not  as  dependent  on  user  input  as  the  inner  loops. 
These  are  described  in  the  following  paragraphs. 

3.6.1.  Electronic  degrees  of  freedom 

The  first  step  in  the  process  of  optimizing  the  EDOF  consists  of 
choosing  a  trial  charge  density  and  a  trial  wavefunction.  In  the  case 
of  the  non-BANDS-type  calculations,  the  trial  wavefunctions  are 
initialized  using  random  numbers,  while  the  trial  charge  density 
is  obtained  from  the  superposition  of  atomic  charge  densities. 
The  BANDS  calculations  are  not  self-consistent,  and  thus  do  not 
feature  a  charge  density  optimization.  In  these  cases  the  charge 
density  obtained  from  the  previously  performed  STATIC  calcula¬ 
tion  is  used  in  the  generation  of  the  starting  wavefunctions. 

Two  iterative  methods  are  used  for  diagonalizing  the  KS 
Hamiltonian:  the  Davidson  blocked  scheme  (DBS)  [57,58],  and 
the  preconditioned  residual  minimization  method-direct  inversion 
in  the  iterative  subspace  (RMM-DIIS)  [10],  Of  the  two,  DBS  is 
known  to  be  the  slower  and  more  stable  option.  Additionally,  the 
subspace  rotation  matrix  is  always  optimized.  These  methods  are 
applied  in  a  manner  that  is  dependent  on  the  calculation  type: 

i.  RELAX  calculations.  Geometry  optimizations  contain  at  least 
one  determination  of  the  system  forces.  The  initial  determi¬ 
nation  consists  of  5  initial  DBS  steps,  followed  by  as  many 
RMM-DHS  steps  as  needed  to  fulfill  the  SEeiec  condition. 
Later  determinations  of  system  forces  are  performed  by  a 
similar  sequence,  but  only  a  single  DBS  step  is  applied  at 
the  outset  of  the  process.  Across  all  databases  the  minimum 
of  number  of  electronic  iterations  for  RELAX  calculations  is 
2.  The  maximum  number  is  set  to  120  for  entries  in  the 
ICSD,  and  60  for  all  others. 

ii.  non-RELAX  calculations.  In  STAT I C  or  BANDS  calculations,  the 
diagonalizations  are  always  performed  using  RMM-DIIS.  The 
minimum  number  of  electronic  iterations  performed  during 
non-RELAX  calculations  is  2,  and  the  maximum  is  120. 

If  the  number  of  iterations  in  the  inner  loop  somehow  exceed 
the  limits  listed  above,  the  calculation  breaks  out  of  this  loop, 
and  the  system  forces  and  energy  are  determined.  If  the  SEion  con¬ 
vergence  condition  is  not  met  the  calculation  re-enters  the  inner 
loop,  and  proceeds  normally. 

Charge  mixing  is  performed  via  Pulay’s  method  [59].  The  imple¬ 
mentation  of  this  charge  mixing  approach  in  the  VASP  package 
depends  on  a  series  of  parameters,  of  which  all  but  the  maximum 
f-quantum  number  handled  by  the  mixer  have  been  left  in  their 
default  state.  This  parameter  is  modified  only  in  systems  included 
in  the  ICSD  database  which  contain  the  elements  listed  in  Tables  4 
and  5.  In  practical  terms,  the  value  applied  in  these  cases  is  the 
maximum  t-quantum  number  found  in  the  PAW  potential,  multi¬ 
plied  by  2. 


3.6.2.  Ionic  degrees  of  freedom  and  lattice  vectors 

The  RELAX  calculation  type  contains  determinations  of  the 
forces  acting  on  the  ions,  as  well  as  the  full  system  stress  tensor. 
The  applied  algorithm  is  the  conjugate  gradients  (CG)  approach 
[60],  which  depends  on  these  quantities  for  the  full  optimization 
of  the  system  geometry,  i.e.  the  ionic  positions,  the  lattice  vectors, 
as  well  as  modifications  of  the  cell  volume.  The  implementation  of 
CG  in  VASP  requires  minimal  user  input,  where  the  only  indepen¬ 
dent  parameter  is  the  initial  scaling  factor  which  is  always  left  at 
its  default  value.  Convergence  of  the  IDOF,  as  stated  above, 
depends  on  the  value  for  the  SEion  parameter,  as  applied  across 
the  various  databases.  The  adopted  Ecut  (see  discussion  on 
“Potentials  and  basis  set”,  Section  3.2)  makes  corrections  for 
Pulay  stresses  unnecessary. 

Forces  acting  on  the  ions  and  stress  tensor  are  subjected  to 
Harris-Foulkes  [61]  corrections.  Molecular  dynamics  based  relax¬ 
ations  are  not  performed  in  the  construction  of  the  databases 
found  in  the  AFLOWLIB  repository,  so  any  related  settings  are  not 
applicable  to  this  work. 

3.7.  Output  options 

The  reproduction  of  the  results  presented  on  the  AFLOWLIB 
website  also  depends  on  a  select  few  parameters  that  govern  the 
output  of  the  DFT  package.  The  density  of  states  plots  are  gener¬ 
ated  from  the  STATIC  calculation.  States  are  plotted  with  a  range 
of  -30  eV  to  45  eV,  and  with  a  resolution  of  5000  points.  The  band 
structures  are  plotted  according  to  the  paths  of  fc-points  generated 
for  a  BANDS  calculation  [18].  All  bands  found  between  -10  eV  and 
10  eV  are  included  in  the  plots. 

4.  Conclusion 

The  AFLOW  standard  described  here  has  been  applied  in  the 
automated  creation  of  the  AFLOWLIB  database  of  material 
properties  in  a  consistent  and  reproducible  manner.  The  use  of 
standardized  parameter  sets  facilitates  the  direct  comparison  of 
properties  between  different  materials,  so  that  specific  trends  can 
be  identified  to  assist  in  the  formulation  of  design  rules  for 
accelerated  materials  development.  Following  this  AFLOW 
standard  should  allow  materials  science  researchers  to  reproduce 
the  results  reported  by  the  AFLOWLIB  consortium,  as  well  as  to 
extend  on  the  database  and  make  meaningful  comparisons  with 
their  own  results. 
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Appendix  A.  Supplementary  material 

Supplementary  data  associated  with  this  article  can  be  found,  in 
the  online  version,  at  http://dx.doi.Org/10.1016/j.commatsci.2015. 
07.019. 
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